# Replication package for 
# "The Role of Foreign Aid in Procuring 
# Civil War Party Consent to Peacekeeping"
# Johannes Karreth, Timothy Passmore, and Jaroslav Tir
# Accepted for publication at Foreign Policy Analysis
# Corresponding author: jkarreth@ursinus.edu

# Please read "0_readme.html" before using this code.

# This file is `analysis.R`
# It generates the main results in the article and in the supporting information.

# The full replication package is posted on Dataverse:
# https://dataverse.harvard.edu/dataverse/jkarreth

## Unit of obs.: conflict

####################################
### Load packages
####################################

rm(list = ls())
library("rio")
library("tidyverse")
library("broom")
library("rstanarm")
library("BayesPostEst")
library("gridExtra")
library("patchwork")
library("texreg")
library("modelsummary")

setwd(dirname(rstudioapi::getSourceEditorContext()$path))


########################
# Load data
########################

unre_1990_ch7_dat <- import("unre_1990_ch7_dat.csv")
unre_1990_noch7_dat <- import("unre_1990_noch7_dat.csv")

########################
# Analyses
########################

# Main analyses -----------------------------------------------------------

# Table 1 Model 1: unrestricted_1990_ch7_gdp_fits_10

unrestricted_1990_ch7_gdp_fits_10 <- stan_glm(unrestricted_binary ~ o_totof_GDP_ihst_last3_lag1_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP + ch7, 
                                              data = unre_1990_ch7_dat,
                                              family = binomial(link = "logit"),
                                              prior = student_t(df = 5, location = 0, scale = 2.5),
                                              prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                              seed = 20230620)

# Table 1 Model 2: unrest_1990_ch7_gdp_fits_10

unrest_1990_ch7_gdp_fits_10 <- stan_glm(unrest_pare ~ o_totof_GDP_ihst_last3_lag1_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP + ch7, 
                                        data = unre_1990_ch7_dat,
                                        family = binomial(link = "logit"),
                                        prior = student_t(df = 5, location = 0, scale = 2.5),
                                        prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                        seed = 20230620)

# write Table 1 to LaTeX:

obs_un_1990_ch7_gdp_reg <- c(dim(unrestricted_1990_ch7_gdp_fits_10$model)[1], dim(unrest_1990_ch7_gdp_fits_10$model)[1])

mcmcReg(mod = list(unrestricted_1990_ch7_gdp_fits_10, unrest_1990_ch7_gdp_fits_10),
        pointest = "median",
        ci = 0.9,
        hpdi = TRUE,
        coefnames = list(c("Intercept", "Total official aid flows/GDP", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points", "Chapter VII mission"),
                         c("Intercept", "Total official aid flows/GDP", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points", "Chapter VII mission")),
        format = "latex",
        custom.model.names = c("Unrestricted vs. any other, including no PKO", "Unrestricted vs. partial/restricted/no consent, PKOs only"),
        digits = 2,
        caption = "Determinants of consent. Columns represent two different comparisons of consent vs. other outcomes. Cases: ceasefires that last at least one month in civil conflicts, 1990-2011, Chapter VII missions included.",
        caption.above = TRUE,
        label = "tab:un_1990_ch7_gdp_reg",
        booktabs = TRUE,
        custom.gof.rows = list(Observations = obs_un_1990_ch7_gdp_reg),
        file = "Output/Tables/table1.tex")

# Table 2 Model 1: unrestricted_1990_noch7_gdp_fits_10

unrestricted_1990_noch7_gdp_fits_10 <- stan_glm(unrestricted_binary ~ o_totof_GDP_ihst_last3_lag1_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP, 
                                                data = unre_1990_noch7_dat,
                                                family = binomial(link = "logit"),
                                                prior = student_t(df = 5, location = 0, scale = 2.5),
                                                prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                                seed = 20230620)

# Table 2 Model 2: unrest_1990_noch7_gdp_fits_10

unrest_1990_noch7_gdp_fits_10 <- stan_glm(unrest_pare ~ o_totof_GDP_ihst_last3_lag1_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP, 
                                          data = unre_1990_noch7_dat,
                                          family = binomial(link = "logit"),
                                          prior = student_t(df = 5, location = 0, scale = 2.5),
                                          prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                          seed = 20230620)

# write Table 2 to LaTeX:

obs_un_1990_noch7_gdp_reg <- c(dim(unrestricted_1990_noch7_gdp_fits_10$model)[1], dim(unrest_1990_noch7_gdp_fits_10$model)[1])

mcmcReg(mod = list(unrestricted_1990_noch7_gdp_fits_10, unrest_1990_noch7_gdp_fits_10),
        pointest = "median",
        ci = 0.9,
        hpdi = TRUE,
        coefnames = list(c("Intercept", "Total official aid flows/GDP", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points"),
                         c("Intercept", "Total official aid flows/GDP", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points")),
        format = "latex",
        custom.model.names = c("Unrestricted vs. any other, including no PKO", "Unrestricted vs. partial/restricted/no consent, PKOs only"),
        digits = 2,
        caption = "Determinants of consent. Columns represent two different comparisons of consent vs. other outcomes. Cases: ceasefires that last at least one month in civil conflicts, 1990-2011, Chapter VII missions excluded.",
        caption.above = TRUE,
        label = "tab:un_1990_noch7_gdp_reg",
        booktabs = TRUE,
        custom.gof.rows = list(Observations = obs_un_1990_noch7_gdp_reg),
        file = "Output/Tables/table2.tex")



# Robustness tests --------------------------------------------------------

# include pre-war aid (Tables A3 and A4)

# T1M1
unrestricted_1990_ch7_gdp_prewar_fits_10 <- stan_glm(unrestricted_binary ~ o_totof_GDP_ihst_last3_lag1_s + o_totof_GDP_ihst_last3_lag1_prewar_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP + ch7, 
                                              data = unre_1990_ch7_dat,
                                              family = binomial(link = "logit"),
                                              prior = student_t(df = 5, location = 0, scale = 2.5),
                                              prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                              seed = 20230620)
# T1M2
unrest_1990_ch7_gdp_prewar_fits_10 <- stan_glm(unrest_pare ~ o_totof_GDP_ihst_last3_lag1_s + o_totof_GDP_ihst_last3_lag1_prewar_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP + ch7, 
                                        data = unre_1990_ch7_dat,
                                        family = binomial(link = "logit"),
                                        prior = student_t(df = 5, location = 0, scale = 2.5),
                                        prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                        seed = 20230620)

# Write pre-war Table 1 to LaTeX: 

obs_un_1990_ch7_gdp_prewar_reg <- c(dim(unrestricted_1990_ch7_gdp_prewar_fits_10$model)[1], dim(unrest_1990_ch7_gdp_prewar_fits_10$model)[1])

mcmcReg(mod = list(unrestricted_1990_ch7_gdp_prewar_fits_10, unrest_1990_ch7_gdp_prewar_fits_10),
        pointest = "median",
        ci = 0.9,
        hpdi = TRUE,
        coefnames = list(c("Intercept", "Total official aid flows/GDP", "Total official aid flows/GDP before war", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points", "Chapter VII mission"),
                         c("Intercept", "Total official aid flows/GDP", "Total official aid flows/GDP before war", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points", "Chapter VII mission")),
        format = "latex",
        custom.model.names = c("Unrestricted vs. any other, including no PKO", "Unrestricted vs. partial/restricted/no consent, PKOs only"),
        digits = 2,
        caption = "Determinants of consent. Columns represent two different comparisons of consent vs. other outcomes. Cases: ceasefires that last at least one month in civil conflicts, 1990-2011, Chapter VII missions included.",
        caption.above = TRUE,
        label = "tab:un_1990_ch7_gdp_prewar_reg",
        booktabs = TRUE,
        custom.gof.rows = list(Observations = obs_un_1990_ch7_gdp_prewar_reg),
        file = "Output/Tables/tableA3_prewar.tex")

# T2M1
unrestricted_1990_noch7_gdp_prewar_fits_10 <- stan_glm(unrestricted_binary ~ o_totof_GDP_ihst_last3_lag1_s + o_totof_GDP_ihst_last3_lag1_prewar_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP, 
                                                data = unre_1990_noch7_dat,
                                                family = binomial(link = "logit"),
                                                prior = student_t(df = 5, location = 0, scale = 2.5),
                                                prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                                seed = 20230620)
# T2M2
unrest_1990_noch7_gdp_prewar_fits_10 <- stan_glm(unrest_pare ~ o_totof_GDP_ihst_last3_lag1_s + o_totof_GDP_ihst_last3_lag1_prewar_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP, 
                                          data = unre_1990_noch7_dat,
                                          family = binomial(link = "logit"),
                                          prior = student_t(df = 5, location = 0, scale = 2.5),
                                          prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                          seed = 20230620)
# write pre-War Table 2 to LaTeX:
obs_un_1990_noch7_gdp_prewar_reg <- c(dim(unrestricted_1990_noch7_gdp_prewar_fits_10$model)[1], dim(unrest_1990_noch7_gdp_prewar_fits_10$model)[1])

mcmcReg(mod = list(unrestricted_1990_noch7_gdp_prewar_fits_10, unrest_1990_noch7_gdp_prewar_fits_10),
        pointest = "median",
        ci = 0.9,
        hpdi = TRUE,
        coefnames = list(c("Intercept", "Total official aid flows/GDP", "Total official aid flows/GDP before war", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points"),
                         c("Intercept", "Total official aid flows/GDP", "Total official aid flows/GDP before war", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points")),
        format = "latex",
        custom.model.names = c("Unrestricted vs. any other, including no PKO", "Unrestricted vs. partial/restricted/no consent, PKOs only"),
        digits = 2,
        caption = "Determinants of consent. Columns represent two different comparisons of consent vs. other outcomes. Cases: ceasefires that last at least one month in civil conflicts, 1990-2011, Chapter VII missions excluded.",
        caption.above = TRUE,
        label = "tab:un_1990_noch7_gdp_prewar_reg",
        booktabs = TRUE,
        custom.gof.rows = list(Observations = obs_un_1990_noch7_gdp_prewar_reg),
        file = "Output/Tables/tableA4_prewar.tex")

# add largest donor as a second explanatory variable (Tables A5 and A6)

# Table 1 Model 1: unrestricted_1990_ch7_gdp_maxdonor_fits_10

unrestricted_1990_ch7_gdp_otherMaxdonor_fits_10 <- stan_glm(unrestricted_binary ~ I(o_totof_GDP_last3_lag1 - o_totof_maxdonor_GDP_last3_lag1) +
                                                              o_totof_maxdonor_GDP_last3_lag1 + 
                                                         gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + 
                                                         major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP + ch7, 
                                                       data = unre_1990_ch7_dat,
                                                       family = binomial(link = "logit"),
                                                       prior = student_t(df = 5, location = 0, scale = 2.5),
                                                       prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                                       seed = 20230620)

# Table 1 Model 2: unrest_1990_ch7_gdp_maxdonor_fits_10

unrest_1990_ch7_gdp_otherMaxdonor_fits_10 <- stan_glm(unrest_pare ~ I(o_totof_GDP_last3_lag1 - o_totof_maxdonor_GDP_last3_lag1) +
                                                        o_totof_maxdonor_GDP_last3_lag1 + 
                                                   gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + 
                                                   major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP + ch7, 
                                                 data = unre_1990_ch7_dat,
                                                 family = binomial(link = "logit"),
                                                 prior = student_t(df = 5, location = 0, scale = 2.5),
                                                 prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                                 seed = 20230620)

# write Table 1 to LaTeX:

obs_un_1990_ch7_gdp_otherMaxdonor_reg <- c(dim(unrestricted_1990_ch7_gdp_otherMaxdonor_fits_10$model)[1], dim(unrest_1990_ch7_gdp_otherMaxdonor_fits_10$model)[1])

mcmcReg(mod = list(unrestricted_1990_ch7_gdp_otherMaxdonor_fits_10, unrest_1990_ch7_gdp_otherMaxdonor_fits_10),
        pointest = "median",
        ci = 0.9,
        hpdi = TRUE,
        coefnames = list(c("Intercept", "Total official aid flows from all sources, excluding largest state donor/GDP", "Total official aid flows from largest state donor/GDP", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points", "Chapter VII mission"),
                         c("Intercept", "Total official aid flows from all sources, excluding largest state donor/GDP", "Total official aid flows from largest state donor/GDP", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points", "Chapter VII mission")),
        format = "latex",
        custom.model.names = c("Unrestricted vs. any other, including no PKO", "Unrestricted vs. partial/restricted/no consent, PKOs only"),
        digits = 2,
        caption = "Determinants of consent. Columns represent two different comparisons of consent vs. other outcomes. Cases: ceasefires that last at least one month in civil conflicts, 1990-2011, Chapter VII missions included.",
        caption.above = TRUE,
        label = "tab:un_1990_ch7_gdp_otherMaxdonor_reg",
        booktabs = TRUE,
        custom.gof.rows = list(Observations = obs_un_1990_ch7_gdp_otherMaxdonor_reg),
        file = "Output/Tables/tableA5_otherMaxdonor.tex")

# Table 2 Model 1: unrestricted_1990_noch7_gdp_maxdonor_fits_10

unrestricted_1990_noch7_gdp_otherMaxdonor_fits_10 <- stan_glm(unrestricted_binary ~ I(o_totof_GDP_last3_lag1 - o_totof_maxdonor_GDP_last3_lag1) +
                                                                o_totof_maxdonor_GDP_last3_lag1 + 
                                                           gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + 
                                                           major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP, 
                                                         data = unre_1990_noch7_dat,
                                                         family = binomial(link = "logit"),
                                                         prior = student_t(df = 5, location = 0, scale = 2.5),
                                                         prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                                         seed = 20230620)

# Table 2 Model 2: unrest_1990_noch7_gdp_maxdonor_fits_10

unrest_1990_noch7_gdp_otherMaxdonor_fits_10 <- stan_glm(unrest_pare ~ I(o_totof_GDP_last3_lag1 - o_totof_maxdonor_GDP_last3_lag1) +
                                                          o_totof_maxdonor_GDP_last3_lag1 + 
                                                     gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + 
                                                     major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP, 
                                                   data = unre_1990_noch7_dat,
                                                   family = binomial(link = "logit"),
                                                   prior = student_t(df = 5, location = 0, scale = 2.5),
                                                   prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                                   seed = 20230620)

# write Table 2 to LaTeX:

obs_un_1990_noch7_gdp_otherMaxdonor_reg <- c(dim(unrestricted_1990_noch7_gdp_otherMaxdonor_fits_10$model)[1], dim(unrest_1990_noch7_gdp_otherMaxdonor_fits_10$model)[1])

mcmcReg(mod = list(unrestricted_1990_noch7_gdp_otherMaxdonor_fits_10, unrest_1990_noch7_gdp_otherMaxdonor_fits_10),
        pointest = "median",
        ci = 0.9,
        hpdi = TRUE,
        coefnames = list(c("Intercept", "Total official aid flows from all sources, excluding largest state donor/GDP", "Total official aid flows from largest state donor/GDP", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points"),
                         c("Intercept", "Total official aid flows from all sources, excluding largest state donor/GDP", "Total official aid flows from largest state donor/GDP", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points")),
        format = "latex",
        custom.model.names = c("Unrestricted vs. any other, including no PKO", "Unrestricted vs. partial/restricted/no consent, PKOs only"),
        digits = 2,
        caption = "Determinants of consent. Columns represent two different comparisons of consent vs. other outcomes. Cases: ceasefires that last at least one month in civil conflicts, 1990-2011, Chapter VII missions excluded.",
        caption.above = TRUE,
        label = "tab:un_1990_noch7_gdp_otherMaxdonor_reg",
        booktabs = TRUE,
        custom.gof.rows = list(Observations = obs_un_1990_noch7_gdp_otherMaxdonor_reg),
        file = "Output/Tables/tableA6_otherMaxdonor.tex")

# specific donors (Tables A7 and A8)

# Table 1 Model 1: unrestricted_1990_ch7_gdp_donors_fits_10

unrestricted_1990_ch7_gdp_donors_fits_10 <- stan_glm(unrestricted_binary ~ o_totof_CAN_GDP_ihst_last3_lag1_s + o_totof_FRN_GDP_ihst_last3_lag1_s + o_totof_GMY_GDP_ihst_last3_lag1_s + o_totof_JPN_GDP_ihst_last3_lag1_s + o_totof_UKG_GDP_ihst_last3_lag1_s + o_totof_USA_GDP_ihst_last3_lag1_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP + ch7, 
                                                     data = unre_1990_ch7_dat,
                                                     family = binomial(link = "logit"),
                                                     prior = student_t(df = 5, location = 0, scale = 2.5),
                                                     prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                                     seed = 20230620)

# Table 1 Model 2: unrest_1990_ch7_gdp_donors_fits_10

unrest_1990_ch7_gdp_donors_fits_10 <- stan_glm(unrest_pare ~ o_totof_CAN_GDP_ihst_last3_lag1_s + o_totof_FRN_GDP_ihst_last3_lag1_s + o_totof_GMY_GDP_ihst_last3_lag1_s + o_totof_JPN_GDP_ihst_last3_lag1_s + o_totof_UKG_GDP_ihst_last3_lag1_s + o_totof_USA_GDP_ihst_last3_lag1_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP + ch7, 
                                               data = unre_1990_ch7_dat,
                                               family = binomial(link = "logit"),
                                               prior = student_t(df = 5, location = 0, scale = 2.5),
                                               prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                               seed = 20230620)

# write Table 1 to LaTeX:

obs_un_1990_ch7_gdp__donors_reg <- c(dim(unrestricted_1990_ch7_gdp_donors_fits_10$model)[1], dim(unrest_1990_ch7_gdp_donors_fits_10$model)[1])

mcmcReg(mod = list(unrestricted_1990_ch7_gdp_donors_fits_10, unrest_1990_ch7_gdp_donors_fits_10),
        pointest = "median",
        ci = 0.9,
        hpdi = TRUE,
        coefnames = list(c("Intercept", "Total official aid flows from Canada/GDP", "Total official aid flows from France/GDP", "Total official aid flows from Germany/GDP", "Total official aid flows from Japan/GDP", "Total official aid flows from UK/GDP", "Total official aid flows from US/GDP", 
                           "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points", "Chapter VII mission"),
                         c("Intercept", "Total official aid flows from Canada/GDP", "Total official aid flows from France/GDP", "Total official aid flows from Germany/GDP", "Total official aid flows from Japan/GDP", "Total official aid flows from UK/GDP", "Total official aid flows from US/GDP", 
                           "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points", "Chapter VII mission")),
        format = "latex",
        custom.model.names = c("Unrestricted vs. any other, including no PKO", "Unrestricted vs. partial/restricted/no consent, PKOs only"),
        digits = 2,
        caption = "Determinants of consent. Columns represent two different comparisons of consent vs. other outcomes. Cases: ceasefires that last at least one month in civil conflicts, 1990-2011, Chapter VII missions included.",
        caption.above = TRUE,
        label = "tab:un_1990_ch7_gdp__donors_reg",
        booktabs = TRUE,
        custom.gof.rows = list(Observations = obs_un_1990_ch7_gdp__donors_reg),
        file = "Output/Tables/tableA7_donors.tex")

# Table 2 Model 1: unrestricted_1990_noch7_gdp_donors_fits_10

unrestricted_1990_noch7_gdp_donors_fits_10 <- stan_glm(unrestricted_binary ~ o_totof_CAN_GDP_ihst_last3_lag1_s + o_totof_FRN_GDP_ihst_last3_lag1_s + o_totof_GMY_GDP_ihst_last3_lag1_s + o_totof_JPN_GDP_ihst_last3_lag1_s + o_totof_UKG_GDP_ihst_last3_lag1_s + o_totof_USA_GDP_ihst_last3_lag1_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP, 
                                                       data = unre_1990_noch7_dat,
                                                       family = binomial(link = "logit"),
                                                       prior = student_t(df = 5, location = 0, scale = 2.5),
                                                       prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                                       seed = 20230620)

# Table 2 Model 2: unrest_1990_noch7_gdp_donors_fits_10

unrest_1990_noch7_gdp_donors_fits_10 <- stan_glm(unrest_pare ~ o_totof_CAN_GDP_ihst_last3_lag1_s + o_totof_FRN_GDP_ihst_last3_lag1_s + o_totof_GMY_GDP_ihst_last3_lag1_s + o_totof_JPN_GDP_ihst_last3_lag1_s + o_totof_UKG_GDP_ihst_last3_lag1_s + o_totof_USA_GDP_ihst_last3_lag1_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP, 
                                                 data = unre_1990_noch7_dat,
                                                 family = binomial(link = "logit"),
                                                 prior = student_t(df = 5, location = 0, scale = 2.5),
                                                 prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                                 seed = 20230620)

# write Table 2 to LaTeX:

obs_un_1990_noch7_gdp__donors_reg <- c(dim(unrestricted_1990_noch7_gdp_donors_fits_10$model)[1], dim(unrest_1990_noch7_gdp_donors_fits_10$model)[1])

mcmcReg(mod = list(unrestricted_1990_noch7_gdp_donors_fits_10, unrest_1990_noch7_gdp_donors_fits_10),
        pointest = "median",
        ci = 0.9,
        hpdi = TRUE,
        coefnames = list(c("Intercept", "Total official aid flows from Canada/GDP", "Total official aid flows from France/GDP", "Total official aid flows from Germany/GDP", "Total official aid flows from Japan/GDP", "Total official aid flows from UK/GDP", "Total official aid flows from US/GDP", 
                           "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points"),
                         c("Intercept", "Total official aid flows from Canada/GDP", "Total official aid flows from France/GDP", "Total official aid flows from Germany/GDP", "Total official aid flows from Japan/GDP", "Total official aid flows from UK/GDP", "Total official aid flows from US/GDP", 
                           "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points")),
        format = "latex",
        custom.model.names = c("Unrestricted vs. any other, including no PKO", "Unrestricted vs. partial/restricted/no consent, PKOs only"),
        digits = 2,
        caption = "Determinants of consent. Columns represent two different comparisons of consent vs. other outcomes. Cases: ceasefires that last at least one month in civil conflicts, 1990-2011, Chapter VII missions excluded.",
        caption.above = TRUE,
        label = "tab:un_1990_noch7_gdp__donors_reg",
        booktabs = TRUE,
        custom.gof.rows = list(Observations = obs_un_1990_noch7_gdp__donors_reg),
        file = "Output/Tables/tableA8_donors.tex")

# Linear probability models (Tables A9 and A10)

# Table 1 Model 1: unrestricted_1990_ch7_gdp_fits_10_lpm

unrestricted_1990_ch7_gdp_fits_10_lpm <- stan_glm(unrestricted_binary ~ o_totof_GDP_ihst_last3_lag1_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP + ch7, 
                                                  data = unre_1990_ch7_dat,
                                                  family = gaussian(),
                                                  prior = student_t(df = 5, location = 0, scale = 2.5),
                                                  prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                                  seed = 20230620)

# Table 1 Model 2: unrest_1990_ch7_gdp_fits_10_lpm

unrest_1990_ch7_gdp_fits_10_lpm <- stan_glm(unrest_pare ~ o_totof_GDP_ihst_last3_lag1_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP + ch7, 
                                            data = unre_1990_ch7_dat,
                                            family = gaussian(),
                                            prior = student_t(df = 5, location = 0, scale = 2.5),
                                            prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                            seed = 20230620)

# write Table 1 to LaTeX:

obs_un_1990_ch7_gdp_reg <- c(dim(unrestricted_1990_ch7_gdp_fits_10_lpm$model)[1], dim(unrest_1990_ch7_gdp_fits_10_lpm$model)[1])

mcmcReg(mod = list(unrestricted_1990_ch7_gdp_fits_10_lpm, unrest_1990_ch7_gdp_fits_10_lpm),
        pointest = "median",
        ci = 0.9,
        hpdi = TRUE,
        pars = "[^sigma]", regex = TRUE,
        coefnames = list(c("Intercept", "Total official aid flows/GDP", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points", "Chapter VII mission"),
                         c("Intercept", "Total official aid flows/GDP", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points", "Chapter VII mission")),
        format = "latex",
        custom.model.names = c("Unrestricted vs. any other, including no PKO", "Unrestricted vs. partial/restricted/no consent, PKOs only"),
        digits = 2,
        caption = "Determinants of consent. Columns represent two different comparisons of consent vs. other outcomes. Cases: ceasefires that last at least one month in civil conflicts, 1990-2011, Chapter VII missions included.",
        caption.above = TRUE,
        label = "tab:un_1990_ch7_gdp_reg_lpm",
        booktabs = TRUE,
        custom.gof.rows = list(Observations = obs_un_1990_ch7_gdp_reg),
        file = "Output/Tables/tableA9_lpm.tex")

# Table 2 Model 1: unrestricted_1990_noch7_gdp_fits_10_lpm

unrestricted_1990_noch7_gdp_fits_10_lpm <- stan_glm(unrestricted_binary ~ o_totof_GDP_ihst_last3_lag1_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP, 
                                                    data = unre_1990_noch7_dat,
                                                    family = gaussian(),
                                                    prior = student_t(df = 5, location = 0, scale = 2.5),
                                                    prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                                    seed = 20230620)

# Table 2 Model 2: unrest_1990_noch7_gdp_fits_10_lpm

unrest_1990_noch7_gdp_fits_10_lpm <- stan_glm(unrest_pare ~ o_totof_GDP_ihst_last3_lag1_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + covIP, 
                                              data = unre_1990_noch7_dat,
                                              family = gaussian(),
                                              prior = student_t(df = 5, location = 0, scale = 2.5),
                                              prior_intercept = student_t(df = 5, location = 0, scale = 2.5), cores = 4,
                                              seed = 20230620)

# write Table 2 to LaTeX:

obs_un_1990_noch7_gdp_reg <- c(dim(unrestricted_1990_noch7_gdp_fits_10_lpm$model)[1], dim(unrest_1990_noch7_gdp_fits_10_lpm$model)[1])

mcmcReg(mod = list(unrestricted_1990_noch7_gdp_fits_10_lpm, unrest_1990_noch7_gdp_fits_10_lpm),
        pointest = "median",
        ci = 0.9,
        hpdi = TRUE,
        pars = "[^sigma]", regex = TRUE,
        coefnames = list(c("Intercept", "Total official aid flows/GDP", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points"),
                         c("Intercept", "Total official aid flows/GDP", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points")),
        format = "latex",
        custom.model.names = c("Unrestricted vs. any other, including no PKO", "Unrestricted vs. partial/restricted/no consent, PKOs only"),
        digits = 2,
        caption = "Determinants of consent. Columns represent two different comparisons of consent vs. other outcomes. Cases: ceasefires that last at least one month in civil conflicts, 1990-2011, Chapter VII missions excluded.",
        caption.above = TRUE,
        label = "tab:un_1990_noch7_gdp_reg_lpm",
        booktabs = TRUE,
        custom.gof.rows = list(Observations = obs_un_1990_noch7_gdp_reg),
        file = "Output/Tables/tableA10_lpm.tex")


# Frequentist logistic regression (Tables A11 and A12)

# Table 1 Model 1: unrestricted_1990_ch7_gdp_fits_10_flogit

unrestricted_1990_ch7_gdp_fits_10_flogit <- glm(unrestricted_binary ~ o_totof_GDP_ihst_last3_lag1_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + ch7, # removed covIP b/c separation
                                                data = unre_1990_ch7_dat,
                                                family = binomial(link = "logit"))

# Table 1 Model 2: unrest_1990_ch7_gdp_fits_10_flogit

unrest_1990_ch7_gdp_fits_10_flogit <- glm(unrest_pare ~ o_totof_GDP_ihst_last3_lag1_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist + ch7, # removed covIP b/c quasi-separation
                                          data = unre_1990_ch7_dat,
                                          family = binomial(link = "logit"))


# write Table 1 to LaTeX:

modelsummary(models = list("Unrestricted vs. any other, including no PKO" = unrestricted_1990_ch7_gdp_fits_10_flogit, "Unrestricted vs. partial/restricted/no consent, PKOs only" = unrest_1990_ch7_gdp_fits_10_flogit),
             vcov = "stata",
             stars = c("*" = .1),
             coef_rename = c("Intercept", "Total official aid flows/GDP", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Chapter VII mission"),
             output = "Output/Tables/tableA11_flogit.tex",
             fmt = 2,
             title = "Determinants of consent. Columns represent two different comparisons of consent vs. other outcomes. Cases: ceasefires that last at least one month in civil conflicts, 1990-2011, Chapter VII missions included.\\label{tab:un_1990_ch7_gdp_reg_flogit}",
             align = "ldd")

# Table 2 Model 1: unrestricted_1990_noch7_gdp_fits_10_flogit

unrestricted_1990_noch7_gdp_fits_10_flogit <- glm(unrestricted_binary ~ o_totof_GDP_ihst_last3_lag1_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist, # removed covIP b/c separation
                                                  data = unre_1990_noch7_dat,
                                                  family = binomial(link = "logit"))

# Table 2 Model 2: unrest_1990_noch7_gdp_fits_10_flogit

unrest_1990_noch7_gdp_fits_10_flogit <- glm(unrest_pare ~ o_totof_GDP_ihst_last3_lag1_s + gdppc_WDI_PW_ln_last3_lag1 + victory + treaty + lndead + lngarm + faction + neighbor + wardur2 + major + polity2_P4_last3_lag1 + gems_drugs_oil + shortestP5dist + longestP5dist, # removed covIP b/c separation
                                            data = unre_1990_noch7_dat,
                                            family = binomial(link = "logit"))

# write Table 2 to LaTeX:

modelsummary(models = list("Unrestricted vs. any other, including no PKO" = unrestricted_1990_noch7_gdp_fits_10_flogit, "Unrestricted vs. partial/restricted/no consent, PKOs only" = unrest_1990_noch7_gdp_fits_10_flogit),
             vcov = "stata",
             stars = c("*" = .1),
             coef_rename = c("Intercept", "Total official aid flows/GDP", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member"),
             output = "Output/Tables/tableA12_flogit.tex",
             fmt = 2,
             title = "Determinants of consent. Columns represent two different comparisons of consent vs. other outcomes. Cases: ceasefires that last at least one month in civil conflicts, 1990-2011, Chapter VII missions excluded.\\label{tab:un_1990_noch7_gdp_reg_flogit}",
             align = "ldd")


# Summary stats & description -----------------------------------------------

models <- list(unrestricted_1990_ch7_gdp_fits_10$model, 
               unrest_1990_ch7_gdp_fits_10$model,
               unrestricted_1990_noch7_gdp_fits_10$model,
               unrest_1990_noch7_gdp_fits_10$model)

sumtab <- list()

for(i in 1:4){
  sumtab[[i]] <- do.call(data.frame, 
                         list(mean = apply(models[[i]], 2, mean, na.rm = TRUE),
                              # median = apply(unrestricted_1990_ch7_gdp_fits_10$model, 2, median, na.rm = TRUE),
                              sd = apply(models[[i]], 2, sd, na.rm = TRUE)#,
                              # min = apply(unrestricted_1990_ch7_gdp_fits_10$model, 2, min, na.rm = TRUE),
                              # max = apply(unrestricted_1990_ch7_gdp_fits_10$model, 2, max, na.rm = TRUE),
                              # n = apply(unrestricted_1990_ch7_gdp_fits_10$model, 2, function(x) length(which(!is.na(x))))
                         ))
}

tab_un_1990_ch7_gdp_sum <- cbind(sumtab[[1]], sumtab[[2]])

tab_un_1990_noch7_gdp_sum <- cbind(sumtab[[3]], sumtab[[4]])

tab_un_1990_ch7_gdp_sum$varname <- c("Consent",
                                     "Total official aid flows/GDP", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points", "Chapter VII mission")

tab_un_1990_noch7_gdp_sum$varname <- c("Consent",
                                       "Total official aid flows/GDP", "GDP per capita (logged)", "Victory", "Settlement", "Total war deaths", "Gov't army size", "More than 2 factions", "Neighbor intervened", "War duration", "Major power involved in war", "Democracy", "Natural resource revenue available", "Distance to closest P5 member", "Distance to farthest P5 member", "Coef of Var on ideal points")

library("xtable")

tab_un_1990_ch7_gdp_sum <- tab_un_1990_ch7_gdp_sum[, c(5, 1:4)]

names(tab_un_1990_ch7_gdp_sum) <- c("Variable", "Mean", "SD", "Mean", "SD")
tab_un_1990_ch7_gdp_sum <- xtable(tab_un_1990_ch7_gdp_sum, col)
caption(tab_un_1990_ch7_gdp_sum) <- "Descriptive statistics for analyses of consent. Analyses use standardized values (original values divided by two standard deviations). Economic variables and democracy scores are measured in the year before the ceasefire that makes conflicts eligible for peacekeeping operations."
print(tab_un_1990_ch7_gdp_sum, 
      file = "Output/Tables/tableA1_sum.tex",
      type = "latex",
      caption.placement = "top",
      size = "scriptsize",
      include.rownames = FALSE,
      booktabs = TRUE)

tab_un_1990_noch7_gdp_sum <- tab_un_1990_noch7_gdp_sum[, c(5, 1:4)]

names(tab_un_1990_noch7_gdp_sum) <- c("Variable", "Mean", "SD", "Mean", "SD")
tab_un_1990_noch7_gdp_sum <- xtable(tab_un_1990_noch7_gdp_sum, col)
caption(tab_un_1990_noch7_gdp_sum) <- "Descriptive statistics for analyses of consent. Analyses use standardized values (original values divided by two standard deviations). Economic variables and democracy scores are measured in the year before the ceasefire that makes conflicts eligible for peacekeeping operations."
print(tab_un_1990_noch7_gdp_sum, 
      file = "Output/Tables/tableA2_sum.tex",
      type = "latex",
      caption.placement = "top",
      size = "scriptsize",
      include.rownames = FALSE,
      booktabs = TRUE)

# Showing variation in aid ---------------------------------

p1 <- ggplot(data = unrestricted_1990_ch7_gdp_fits_10$model, aes(x = o_totof_GDP_ihst_last3_lag1_s)) + 
  geom_histogram(fill = "white", color = "black") + 
  xlab("Total official aid flows/GDP (standardized)") + ylab("Observations in estimation sample") + 
  labs(title = "Sample from Table 1, Model 1") +
  theme_minimal()

p2 <- ggplot(data = unrest_1990_ch7_gdp_fits_10$model, aes(x = o_totof_GDP_ihst_last3_lag1_s)) + 
  geom_histogram(fill = "white", color = "black") + 
  xlab("Total official aid flows/GDP (standardized)") + ylab("Observations in estimation sample") + 
  labs(title = "Sample from Table 1, Model 2") +
  theme_minimal()

p3 <- ggplot(data = unrestricted_1990_noch7_gdp_fits_10$model, aes(x = o_totof_GDP_ihst_last3_lag1_s)) + 
  geom_histogram(fill = "white", color = "black") + 
  xlab("Total official aid flows/GDP (standardized)") + ylab("Observations in estimation sample") + 
  labs(title = "Sample from Table 2, Model 1") +
  theme_minimal()

p4 <- ggplot(data = unrest_1990_noch7_gdp_fits_10$model, aes(x = o_totof_GDP_ihst_last3_lag1_s)) + 
  geom_histogram(fill = "white", color = "black") + 
  xlab("Total official aid flows/GDP (standardized)") + ylab("Observations in estimation sample") + 
  labs(title = "Sample from Table 2, Model 2") +
  theme_minimal()

library("patchwork")
ps <- p1+p2+p3+p4
ggsave(ps, filename = "Output/Graphs/figA1_aid_histograms.pdf", width = 8, height = 8)


########################
# Save all results
########################

save.image(file = "final_workspace.RData")
